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Cosmic rays are a fundamental source of ionization for molecular and diffuse clouds, influencing 
their chemical, thermal, and dynamical evolution. The amount of cosmic rays inside a cloud 
also determines the 7 -ray flux produced by hadronic collisions between cosmic rays and cloud 
material. We study the spectrum of cosmic rays inside and outside of a diffuse cloud, by solving 
the stationary transport equation for cosmic rays including diffusion, advection and energy losses 
due to ionization of neutral hydrogen atoms. We found that the cosmic ray spectrum inside a 
diffuse cloud differs from the one in the interstellar medium (ISM) for energies smaller than 
Ej,,- Ri 100 MeV, irrespective of the model details. Below Eh^, the spectrum is harder (softer) 
than that in the ISM if the latter is a power law with s larger (smaller) than ^ 0.42. As 

a consequence also the ionization rate due to CRs is strongly affected. Assuming an average 
Galactic spectrum similar to the one inferred from AMS-2 and Voyager 1 data, we discuss the 
resulting ionization rate in a typical diffuse cloud. 
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1. Introduction 

The amount of penetration of cosmic rays (CRs) into molecular clouds (MCs) regulates the 
ionization level of clouds and dense cores (for reviews see [||, ^ ||]) and thus affects their dynamical 
evolution and the process of star formation. Moreover, the exclusion of CRs from MCs can reduce 
their gamma-ray emission [Q, js]], which results from the decay of neutral pions produced in inelastic 
interactions of CRs in the dense gas (see for a review). Therefore, it is of prime importance to 
understand wether CRs do penetrate or not MCs. 

The difficulty of modeling the problem of CR penetration into clouds resides in its highly non¬ 
linear nature: CRs generate magnetic turbulence at the cloud border due to streaming instability. 
The level of the turbulence in turn determines the diffusion coefficient of CRs and thus, presumably, 
their capability of penetrate clouds. Some discrepancy exists in the literature, different theoretical 
approaches to the problem giving different results. According to early papers, CRs with energies 
below tens or hundreds of MeV are effectively excluded from MCs [^, while in a more recent 
work, 0 found out that the CR intensity is only slightly reduced inside clouds. However, a direct 
comparison of the two approaches is not straightforward, since the former are kinetic approaches, 
while the latter a two fluid ones. An implicit assumption in all these papers is the fact that streaming 
instability would enhance the magnetic turbulence and cause the exclusion of CRs from clouds. 

Here, we present a solution of the steady-state kinetic transport equation of CRs along a mag¬ 
netic flux tube that encompasses a MC. We generalize the simplified two-zones (in- and out-side 
of the cloud) kinetic approaches by [Q] and by considering the full spatially dependent equa¬ 
tion. Remarkably, we find that the exclusion of CRs from diffuse clouds of typical column density 
Nh ~ 3 X 10^^ cm^^ is effective below an energy of « 100 MeV, independently on the presence 
or not of streaming instability. In fact, the exclusion energy Ehr depends only (and quite weakly) 
on the physical parameters that characterize the ISM and the gas in the cloud. This result suggests 
that: i) the suppression of the gamma-ray emission from a cloud due to CR exclusion is not signif¬ 
icant (the threshold for neutral pion production being Ri 280 MeV), and ii) the intensity of CRs is 
suppressed inside MCs at the particle energies which are most relevant for ionization. 

2. The model 

Consider a cloud of size Lc and of hydrogen density uh threaded in a magnetic field of intensity 
Bq, oriented along the x-axis. Such a one-dimensional configuration is realistic if one considers 
spatial scales smaller than the magnetic field coherence length in the ISM, i.e. ~ 50-100 pc. 
Moreover, the magnetic field strength is assumed to be spatially constant and not to change across 
the transition from outside to inside of the cloud. This is supported by observations [||] showing 
that the magnetic field strength is independent on the density of the interstellar medium as long as 
the latter remains smaller than « 300 cm^^. The cloud is assumed to be immersed in a diffuse, 
hot, and fully ionized ISM of density Following [§- [01> and [^], CRs are assumed to propagate 
along the magnetic field lines only, i.e. diffusion perpendicular to field lines is set to zero. 

For simplicity, the transition between the low density and ionized ISM and the dense cloud is 
taken to be sharp and located at x = (see Figure |T]). In order to determine the CR profile along the 
X— axis we solve the one dimensional steady state equation for the transport of CRs in presence of 


2 



Cosmic ray penetration in dijfuse clouds 


Giovanni Morlino 



1 

1 

1 

1 

Lc 

1 

1 

1 

1 


zone 1 

1 

' zone 2 

1 

1 

zone 3 
Cloud 

1 

zone 2 ' 

1 

1 

zone 1 


1 


1 



1 

1 

1 


1 

1 

1 

Bo 

1 1 

1 1 

1 1 


—OO 0 Xc Xc+Lc 00 


Figure 1: Sketch of the simplified 1-D model used to describe the cloud geometry. 


diffusion, advection and energy losses. Written in the rest frame of the plasma the equation reads: 
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( 2 . 1 ) 


where f{x,p) is the CR distribution function, D{x,p) is the diffusion coefficient, p{x,p) is the 
momentum losses and va {x) = R(;c)^/\/47rp,- is the Alfven speed in the regular magnetic field B{x). 
Nolice fhal in fhis expression we include only fhe ions mass densify, p,-, because al fhe wavelengfhs 
relevanf in fhis work, fhe wave frequency is smaller fhan fhe charge-exchange frequency befween 
ions and neufrals, hence, while ions oscillate with waves, neutrals have no time to couple with 
them (e.g. 1W]). Note that using Eq. we assume a diffusive behavior to CRs both outside and 
inside of the cloud. The validity of such approximation inside the cloud is questionable, because 
the magnetic turbulence is damped very efficiently by the ion-neutral friction. Nevertheless for the 
moment we assume that the propagation is diffusive also inside the cloud. We will show that our 
prediction on the CR spectrum is not strongly affected by this assumption. 

Energy losses are presents only inside the cloud (region 3 in Eig.|i|). In the energy range 
between 100 keV and 1 GeV the loss time T/ due to ionization of neutral hydrogen can be well 
approximated by a power law in momentum, which reads: 


= - = To 
P 


0.1 nipC 


riH 


cm 


-3 


( 2 . 2 ) 


where the normalization and the slope are Tq = 1-46 • 10^ yr and a = 2.58, respectively, and have 
been obtained fitting the energy losses provided by [ [II| ] (see their Eigure 7). 

Eet us now discuss the boundary conditions. Ear from the cloud {x ^ Xc) we impose that 
the CR distribution reduces to the Galactic one fo{p), while at the cloud centre we impose the 
symmetry condition dxf{x,p)\x=Xc+c /2 = 0- The latter condition is different from what has been 
done in [^], where a condition on the CR gradient far from the cloud, rather than at the cloud center 
was imposed, to match the CR gradient that one would expect in a Galactic environment. However, 
such a choice may lead to an unphysical solution, because it is what happens inside the cloud 
that determines the value of the gradient outside of the cloud, and not viceversa. Eor this reason 
we decided to use the symmetry condition, which is valid as long as the diffusion approximation 
holds. In our approach the value of dxf away from the cloud is an outcome of the calculation. We 
also notice that the symmetry condition is valid only for isolated clouds and breaks if the cloud 
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is located near a CR source. Such a situation requires different boundary conditions and will be 
considered elsewhere. 


If we introduce the function g = Ddxf, Eq.( |2.1| ) reduces to a linear differential equation of 
the first order: 

dxg-gVA/D + Q = 0, (2.3) 

The nonlinearity of the problem has been hidden in the function Q{x,p) which plays the role of a 
source/sink term: 

n fivA a f ^ r- , 

(2.4) 
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The solution for g is: 
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D{y,p) 

We can now write the solution above for the simplified geomefry of fhe MC skefched in Fig¬ 
ure ||. Inside fhe MC we assume a consfanf densify of neufral hydrogen of tin = 100 cm^^, wifh a 


ionization fraction of 10^^. For the diffuse ISM we take = 10^^ cm^-*, so that the Alfven speed is 
constant across the transition between the diffuse medium and the MC. For the moment we ignore 
the effect of streaming instability and we assume a Kolmogorov diffusion coefficient outside of the 
MC: D{x,p) = Dkoi{p) ~ I0^^(p/mc)^/^j8 cm^/s, with jS = Vp/c. Inside the MC Alfven waves are 
heavily damped due to ion-neutral friction and the CR diffusion coefficient is Dc S> Dkoi- 

Under these assumptions it is straightforward to derive from Fq. (^) an expression for / 
outside of the cloud: 


f{x,p)= h{p) + —e 
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The exponent outside the integral tells that the CR density outside of the cloud is affected up to a 
distance of: 
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(2.7) 


which can be much larger than the cloud size, and, for a particle energy of « 100 MeV is of the 
order of 100 pc, comparable to the magnetic field’s coherence length in the interstellar medium. 
For distances to the MC larger than Xc the CR distribution reduces to the Galactic one. Strictly 
speaking Fq. ( |2.7[ ) implies that, for F > 100 MeV the 1-D approximation breaks down and a more 
complex transport model should be adopted. However, as we demonstrated in [p^], when the effect 
of the streaming instability is taken into account the validity of the 1-D approach is guaranteed up 
to particle energies of few hundreds of MeV. 

Fq. can be further simplified by assuming fhaf fhe diffusion coefficienf inside fhe cloud is 
Dc ^ LcVa (a condifion which is easily fulfilled). This implies fhaf fhe exponenf inside fhe integral 
reduces fo a 5-function. We also note fhaf a critical momenfum exisfs given by fhe condition 
fhaf fhe loss time is longer than the CR propagation time across the cloud T/(p*) > Ll/2Dc{p*). 
Then, for p^ p*, the CR spatial distribution inside the MC is roughly constant and after some 
manipulations Fq. (2.6) becomes: 

f{x,p) =Mp)- 
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Figure 2: Left panel: spatial profile of the CR density outside of the MC according to the solution (2.14) 
(solid line) for p = pbr ~ 100 MeV /c. The dot-dashed line accounts for the presence of streaming instability 
(see [Q). Right panel: Spectra of CRs inside the cloud (solid lines) assuming that the spectra far for 
the cloud are given by power law in momentum of slope s (dashed lines). The dot-dashed line shows the 
eigenfunction of the problem. 


where we also made use of the continuity of the CR distribution function at the MC border fc = 

We note that Eq. (2.8) can be used to find the solution of the problem also when free streaming 
of CRs (instead of diffusion) is assumed inside the MC. This is because also under such assump¬ 
tion, a momentum p* exists above which the CR distribution function can be considered spatially 
constant. In this case p* is determined by T;(p*) = Lc/vst where v.^ ~ Vp/3 is the free streaming 
velocity of CRs inside the MC. For Lc = 10 pc and uh = 100 cm^^ we get E* =E(p*) « 0.6 MeV. 
The only difference with respect to the diffusive case is that under the assumption of free streaming 
inside the MC the continuity condition f{x^) = /(x+) has to be considered only as an approximate 
one. 


After using Eq. ( |2.2| ) and performing the derivative in momentum, Eq. ( |2.8[ ) becomes: 

f{x,p) = Mp) + [(3 “ (x)fcEpdpfc\e^ ■ 


vaT{p) 


(2.9) 


and can be solved for x = to give: 
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where we used the fact that t is a power law. The adimentional function T] is defined as 

The condition q (p) = 1 defines a breaking momentum pbr> or a breaking energy Ebr which reads 
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In fact, for p S> Pbr (i-S- V{p) ^ 1)> the exponential function in Eq. (2.10) reduces to si 

p/ri{p)5{p — pby). Asa consequence fdp) ~ fo{p)- On the other hand for p <C Pbr (i-C- t](p) <C 
1), the exponential function can be approximated as ~ Q(^p — this latter case, if 

we assume a power law for the Galactic CR spectrum, /o(p) p^\ we get 


fc{p) 


T1(P) 

(s-3)p^ 


[plr '-P^ 1 for P ^ Pbr- 


(2.13) 


which implies that, for p <C Pbr the solution is a power law fc{p) for 5 < 3 and fdp) 

for s>3. Remarkably, when the CR spectrum in the ISM is a power-law in momentum, fo{p) p^^ 
with s = a — 3 ^ 0.42, the slope of the spectrum of CR is identical inside and outside of the cloud. 
In other words the function /o oc p^°-42 i^ eigenfunction of the problem, giving fc=fo- The right 


panel of Figure ^ shows few examples for fdp) obtained from Eq. (2.10) using a simple power law 


for fo°^ p ^ (showed with dashed lines). 


A remarkable property of solution (2.10) is that fc does not depend on the diffusion coefficient 
but only on the Alfven speed in the ISM and on the cloud properties. On the other hand, it is easy 
to show that the CR spectrum outside the MC does depend on the diffusion coefficient as: 


f{x,p) = Mp) + [fdp) - Mp)] . 


(2.14) 


which is shown in Figure ^ (left panel) with a solid line, for a spectrum of Galactic CRs /o(p) 


- 4.7 


p "■' normalized to an energy density of 1 eV cm •*. The typical length scale where the spectrum 


is affected is Xc as defined in Eq.(^^. As demonsfrafed in [1^] when fhe sfreaming insfabilify is 
faken info accounf, fhe solution inside fhe cloud, fc, remains unalfered, while fhe value of my be 
reduced up fo one order of magnifude, depending on fhe level of damping of magnetic furbulence 
presenf in fhe medium oufside fhe cloud. The leff panel of Figure ^ also reporfs fhe solufion when 
sfreaming insfabilify is faken info accounf, buf wifhouf any damping. 


3. Ionization rate induced by CR 


We showed in fhaf, as a consequence of propagafion, fhe CR specfrum inside a MC is 
affecfed for energies below E^r given by Eq.(2.12). This fad has very imporfanf consequences on 
fhe ionization level produced by CRs inside a MC. In fad fhe ionization cross section peaks for 
an incidenf profon energy ~ O.I MeV, and drops as for E > Fpeak- Hence any process able fo 
modify fhe CR specfrum in fhe MeV region will sfrongly affecl fhe ionization degree of chemicals 
inside fhe cloud. The ionization rale is defined as: 

/'£^max 


^ =4;r 


jdEk)or{Ek)dEk, 


(3.1) 


where jk is fhe CR flux of specie k, is fhe ionization cross seclion and I is fhe ionizalion 
potential of fhe considered chemical specie. Here we limif our considerafions fo fhe ionizalion 
rate produced by CR protons on fhe H 2 molecule, , in a typical diffuse cloud. We follow fhe 
calculalion as presenled in [|^] which also include fhe conlribulion fo ionizalion due fo secondary 
eledrons. The Galactic CR specfrum fo use in Eq. ( ^71] ) is nof well know, mainly because fhe local 
measured CR specfrum is affecled by solar modulation. For Ihis reason, following [p^], we adopl 
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Figure 3: Shown are the maximum and minimum Galactic CR proton spectra (solid thick and thin lines, re¬ 
spectively) and the corresponding spectra inside a diffuse cloud with a column density 3 x 10^^ cm^ (dashed 
thick and thin lines). The vertical dot dashed line mark the critical and the breaking energies E* = 0.58 MeV 
and Zibr = 99 MeV. 


a minimum and a maximum spectra which reasonably bound the Galactic spectrum. The lower 
bound is obtained fitting the data collected by the Voyager 1 spacecraft from a region beyond the 
Solar termination shock, in the energy range between ~ 5 MeV and ~ 50 MeV. The CR proton flux 
can be parametrized as 

jp{E) = 2.4 X 10^^{E + Eoy^ eV^'cm^^s^isj.-! ^ (•3 2) 


where the minimum spectrum has 5 = 0.3 and j3 =3.0, while the maximum spectrum has 5 = —0.5 
and j3 = 2.2. The crossover energy is Eq = 500 MeV; at £ S> Eq, both the spectra have a power 
low behavior jp{E) oc E^^^ = E^^-^, as detected by the AMS 02 experiment. For the sake of com¬ 
pleteness we remember that the flux in energy as reported in Eq.(^^ is related to the distribution 
function in momentum, f{p), by the relation: j{E) = fE{E)v{E)/{4n) = v{E)pEf{p)/c^, where 
fsiE) is the distribution function in energy. Using Eq.(|3.2D we can estimate the spectrum, /c, in¬ 
side a cloud. The results are reported in Eigure ^ where we compare the energy spectra fE.mm{E) 
and /£,max(£') with fc,E,mm{E) and fc,E,max{E), assuming = 10 pc and tin = 100 cm^^ which 
implies a column density of Nh — 3 • 10^' cm^. One can see from the figure that already at Ebi = 99 
MeV both the spectra are attenuated by a factor ~ 5, while at E* = 0.58 MeV, the attenuation 
reaches roughly 3 orders of magnitudes for both models. Eor energies below E* our simplified 
solution (2.10D provides incorrect prediction. The ionization rate calculated using the unmodified 


maximum and minimum spectra are 3.59 x 10 


-16 ,-l 


and 3.49 X 10 s respectively. When 


we account for the propagation, the resulting ionization rate drops to 2.61 x 10 
maximum model and to 1.05 x 10^'^ s^', for minimum one. 


-17 


for the 


It is worth noticing that for typical diffuse clouds with column densities ranging between 10^' 
and 10^^ cm^, the measured value of is between 10^'^ and 10^^^ s^' (see, e.g. [|rT|]). This 
implies that CR protons alone cannot account for the whole ionization rate and that the role of 
electrons has to be carefully considered. 
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4. Discussion and conclusions 

In this paper we demonstrated that, due to the balance between advective flux of CRs into the 
MC and energy losses into the cloud, an equilibrium spectrum forms inside the MC, characterized 
by a feature at an energy of the order of Ehr ~ 100 MeV for diffuse MCs. Below E},r the spectrum 
falls below (rises above) the CR spectrum in the ISM if the latter is steeper (harder) than oc 
This fact will have an impact on the estimates of the CR ionization rates in MCs (see 0 for a 
review) which are often computed by assuming a ballistic propagation of CRs with a boundary 
condition fc = fo (i.e. the CR spectrum at the border of the cloud is equal to that of interstellar 
CRs). However, before drawing firm conclusions on this issue, the role of CR electrons must be 
assessed, and this will be done in a forthcoming publication. 

The results obtained here are based on the assumption of stationarity. This assumption is 
justified because the typical time scale of the problem can be estimated from dimensional analysis 
as D/v^ which for a particle energy of 100 MeV gives « 10^(D/D^o/)(vA/100km/s)^^ yr, which is 
always shorter than the dynamical (free-fall) time of the cloud (Gp)^^/^ 10^(n///100 cm^^) yr. 

Finally, we note that a break in the spectrum in the 100 MeV range would not affect signifi¬ 
cantly the gamma-ray luminosity of the cloud, because the threshold for neutral pion production in 
proton-proton interactions is at a larger energy, namely « 280 MeV. 
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